package fem2d;

import math2.FunctionRnToR;
import fem2.AbstractStructuralStaticDemo;
import fem2.Constraint;
import fem2.Debugger;
import fem2.MaterialModel;
import fem2.Mesh;
import fem2.MeshGenerator;
import fem2.MeshUtilities;
import fem2.Model;
import fem2.analysis.Analysis;
import fem2.analysis.NonlinearStaticAnalysis;
import fem2.enu.EchoLevelType;
import fem2.enu.State;
import fem2.material.ExponentialSofteningRule;
import fem2.material.FluencyCriteria;
import fem2.material.IsotropicDamageMMimplexNonlocal;
import fem2.mesh_and_topology.IntegrationPointIndexer;
import fem2.observer.TimeStepObserver;
import fem2.pre_and_post.GidMeshGenerator;
import fem2.strategies.Strategy;

public class LshapeIsotropicDamageImplexNonlocal extends LshapeIsotropicDamage {

	public LshapeIsotropicDamageImplexNonlocal() {
		projectDir = "/home/hbui/kratos_janosch";
		projectName = "LshapeIsotropicDamage.gid";
		// meshFileName = "structured_mesh-implex-nonlocal.msh";
		meshFileName = "q4mesh2-implex-nonlocal.msh";
	}

	@Override
	public Mesh createMesh() {
		String fn = projectDir + '/' + projectName + '/' + meshFileName;
		MeshGenerator mg = new GidMeshGenerator(fn);
		Mesh mesh = mg.getMesh(2);
		return mesh;
	}

	@Override
	public Model createConditions(Model m) {
		Mesh mesh = m.getMesh();
		Constraint c1 = new Constraint(false, false);
		Constraint c2 = new Constraint(true, false);
		// c2.setValue(1, -0.2);
		// c2.setValue(1, -1.0);
		c2.setValue(1, -10.0);
		mesh.setConstraint(c1, MeshUtilities.seekNodesOnSurface(mesh, 0, 1, -1000));
		MeshUtilities.seekNode(mesh, 1000, 0).setConstraint(c2);
		return m;
	}

	@Override
	public MaterialModel createMaterial(Model m) {
		double E = 25850;
		double nu = 0.18;
		double t = 200;
		double ft = 2.7;
		// double Gf = 0.065;
		double e0 = ft / Math.sqrt(E);
		// double ef = Gf / ft;
		double ef = 1e2 * e0;
		FluencyCriteria softeningRule = new ExponentialSofteningRule(e0, ef);

		final double interactionRadius = 50;
		FunctionRnToR weightFunction = new FunctionRnToR() {

			@Override
			public double valueAt(double... x) {
				if (x[0] < interactionRadius && x[0] >= 0.0) {
					return Math.pow(1 - Math.pow(x[0] / interactionRadius, 2), 2);
				} else {
					return 0;
				}
			}

			@Override
			public double[] gradientAt(double... x) {
				return null;
			}
		};

		return new IsotropicDamageMMimplexNonlocal(E, nu, 0, t, State.PLANE_STRAIN, softeningRule,
				e0, weightFunction, interactionRadius, new IntegrationPointIndexer(m.getMesh()));

	}

	@Override
	public void run(int nt) {
		Model m = createModel();

		Analysis an = new NonlinearStaticAnalysis(m, nt);
		Strategy s = an.getStrategy();
		s.setEchoLevel(EchoLevelType.OUTPUT, EchoLevelType.AUX1);
		// ((NewtonRaphsonStrategy) s)
		// .setLineSearchCriteria(LineSearchCriteria.INTERPOLATED_LINE_SEARCH);

		addObservers(m, s, new TimeStepObserver(s));

		an.run();

		// postProcess(m);

		Debugger.warn("Analysis completed");
	}

	public static void main(String[] args) {
		AbstractStructuralStaticDemo demo = new LshapeIsotropicDamageImplexNonlocal();

		// demo.run(10); //for du = -0.2
		// demo.run(40);
		demo.run(200);
		// demo.run(1000);

	}
}
